home *** CD-ROM | disk | FTP | other *** search
/ Tech Arsenal 1 / Tech Arsenal (Arsenal Computer).ISO / tek-19 / iritsm3s.zip / CBZR_AUX.C < prev    next >
C/C++ Source or Header  |  1991-05-19  |  10KB  |  280 lines
  1. /******************************************************************************
  2. * CBzr-Aux.c - Bezier curve auxilary routines.                      *
  3. *******************************************************************************
  4. * Written by Gershon Elber, Mar. 90.                          *
  5. ******************************************************************************/
  6.  
  7. #ifdef __MSDOS__
  8. #include <stdlib.h>
  9. #endif /* __MSDOS__ */
  10.  
  11. #include <ctype.h>
  12. #include <stdio.h>
  13. #include <string.h>
  14. #include "cagd_loc.h"
  15.  
  16. /******************************************************************************
  17. * Given a bezier curve - subdivide it into two at the given parametric value. *
  18. * Returns pointer to first curve in a list of two curves (subdivided ones).   *
  19. * The subdivision is exact result of evaluating the curve at t using the      *
  20. * recursive algorithm - the left resulting points are the left curve, and the *
  21. * right resulting points are the right curve (left is below t).              *
  22. ******************************************************************************/
  23. CagdCrvStruct *BzrCrvSubdivAtParam(CagdCrvStruct *Crv, CagdRType t)
  24. {
  25.     CagdBType
  26.     IsNotRational = !CAGD_IS_RATIONAL_CRV(Crv);
  27.     int i, j, l,
  28.     k = Crv -> Length,
  29.     MaxCoord = CAGD_NUM_OF_PT_COORD(Crv -> PType);
  30.     CagdRType
  31.     t1 = 1.0 - t;
  32.     CagdCrvStruct
  33.     *RCrv = BzrCrvNew(k, Crv -> PType),
  34.     *LCrv = BzrCrvNew(k, Crv -> PType);
  35.  
  36.     /* Copy Curve into RCrv, so we can apply the recursive algo. to it.      */
  37.     for (i = 0; i < k; i++)
  38.     for (j = IsNotRational; j <= MaxCoord; j++)
  39.         RCrv -> Points[j][i] = Crv -> Points[j][i];
  40.  
  41.     for (j = IsNotRational; j <= MaxCoord; j++)
  42.     LCrv -> Points[j][0] = Crv -> Points[j][0];
  43.  
  44.     /* Apply the recursive algorithm to RCrv, and update LCrv with the         */
  45.     /* temporary results. Note we updated the first point of LCrv above.     */
  46.     for (i = 1; i < k; i++) {
  47.     for (l = 0; l < k - i; l++)
  48.         for (j = IsNotRational; j <= MaxCoord; j++)
  49.         RCrv -> Points[j][l] = RCrv -> Points[j][l] * t1 +
  50.                      RCrv -> Points[j][l+1] * t;
  51.     /* Copy temporary result to LCrv: */
  52.     for (j = IsNotRational; j <= MaxCoord; j++)
  53.         LCrv -> Points[j][i] = RCrv -> Points[j][0];
  54.     }
  55.     LCrv -> Pnext = RCrv;
  56.  
  57.     return LCrv;
  58. }
  59.  
  60. /******************************************************************************
  61. * Return a new curve, identical to the original but with one degree higher    *
  62. * Let old control polygon be P(i), i = 0 to k-1, and Q(i) be new one then:    *
  63. *               i        k-i                          *
  64. * Q(0) = P(0), Q(i) = --- P(i-1) + (---) P(i), Q(k) = P(k-1).              *
  65. *               k         k                          *
  66. ******************************************************************************/
  67. CagdCrvStruct *BzrCrvDegreeRaise(CagdCrvStruct *Crv)
  68. {
  69.     CagdBType
  70.     IsNotRational = !CAGD_IS_RATIONAL_CRV(Crv);
  71.     int i, j,
  72.     k = Crv -> Length,
  73.     MaxCoord = CAGD_NUM_OF_PT_COORD(Crv -> PType);
  74.     CagdCrvStruct
  75.     *RaisedCrv = BzrCrvNew(k + 1, Crv -> PType);
  76.  
  77.     for (j = IsNotRational; j <= MaxCoord; j++)                /* Q(0). */
  78.     RaisedCrv -> Points[j][0] = Crv -> Points[j][0];
  79.  
  80.     for (i = 1; i < k; i++)                        /* Q(i). */
  81.     for (j = IsNotRational; j <= MaxCoord; j++)
  82.         RaisedCrv -> Points[j][i] =
  83.         Crv -> Points[j][i-1] * (i / ((CagdRType) k)) +
  84.         Crv -> Points[j][i] * ((k - i) / ((CagdRType) k));
  85.  
  86.     for (j = IsNotRational; j <= MaxCoord; j++)                /* Q(k). */
  87.     RaisedCrv -> Points[j][k] = Crv -> Points[j][k-1];
  88.  
  89.     return RaisedCrv;
  90. }
  91.  
  92. /******************************************************************************
  93. * Return a normalized vector, equal to the tangent to Crv at parameter t.     *
  94. * Algorithm: pseudo subdivide Crv at t and using control point of subdivided  *
  95. * curve find the tangent as the difference of the 2 end points.              *
  96. ******************************************************************************/
  97. CagdVecStruct *BzrCrvTangent(CagdCrvStruct *Crv, CagdRType t)
  98. {
  99.     static CagdVecStruct P2;
  100.     CagdVecStruct P1;
  101.     CagdBType
  102.     IsNotRational = !CAGD_IS_RATIONAL_CRV(Crv);
  103.     int i, j, l,
  104.     k = Crv -> Length,
  105.     MaxCoord = CAGD_NUM_OF_PT_COORD(Crv -> PType);
  106.     CagdRType
  107.     t1 = 1.0 - t;
  108.     CagdCrvStruct *RCrv;
  109.  
  110.     if (APX_EQ(t, 0.0)) {
  111.     /* Use Crv starting tangent direction. */
  112.     CagdCoerceToE3(P1.Vec, Crv -> Points, 0, Crv -> PType);
  113.     CagdCoerceToE3(P2.Vec, Crv -> Points, 1, Crv -> PType);
  114.     }
  115.     else if (APX_EQ(t, 1.0)) {
  116.     /* Use Crv ending tangent direction. */
  117.     CagdCoerceToE3(P1.Vec, Crv -> Points, k - 2, Crv -> PType);
  118.     CagdCoerceToE3(P2.Vec, Crv -> Points, k - 1, Crv -> PType);
  119.     }
  120.     else {
  121.     RCrv = BzrCrvNew(k, Crv -> PType);
  122.  
  123.     /* Copy Crv into RCrv, so we can apply the recursive algo. to it. */
  124.     for (i = 0; i < k; i++)
  125.         for (j = IsNotRational; j <= MaxCoord; j++)
  126.         RCrv -> Points[j][i] = Crv -> Points[j][i];
  127.  
  128.     /* Apply the recursive algorithm to RCrv. */
  129.     for (i = 1; i < k; i++)
  130.         for (l = 0; l < k - i; l++)
  131.         for (j = IsNotRational; j <= MaxCoord; j++)
  132.             RCrv -> Points[j][l] = RCrv -> Points[j][l] * t1 +
  133.                        RCrv -> Points[j][l+1] * t;
  134.  
  135.     CagdCoerceToE3(P1.Vec, RCrv -> Points, 0, RCrv -> PType);
  136.     CagdCoerceToE3(P2.Vec, RCrv -> Points, 1, RCrv -> PType);
  137.  
  138.     CagdCrvFree(RCrv);
  139.     }
  140.  
  141.     CAGD_SUB_VECTOR(P2, P1);
  142.  
  143.     CAGD_NORMALIZE_VECTOR(P2);                /* Normalize the vector. */
  144.  
  145.     return &P2;
  146. }
  147.  
  148. /******************************************************************************
  149. * Return a normalized vector, equal to the binormal to Crv at parameter t.    *
  150. * Algorithm: pseudo subdivide Crv at t and using 3 consecutive control point  *
  151. * (p1, p2, p3) of subdivided curve find the bi-normal as the cross product of *
  152. * the difference (p1 - p2) and (p2 - p3).                      *
  153. *   Since a curve may have not BiNormal at inflection points or if the 3      *
  154. * points are colinear, NULL will be returned at such cases.              *
  155. ******************************************************************************/
  156. CagdVecStruct *BzrCrvBiNormal(CagdCrvStruct *Crv, CagdRType t)
  157. {
  158.     static CagdVecStruct P3;
  159.     CagdVecStruct P1, P2;
  160.     CagdBType
  161.     IsNotRational = !CAGD_IS_RATIONAL_CRV(Crv);
  162.     int i, j, l,
  163.     k = Crv -> Length,
  164.     MaxCoord = CAGD_NUM_OF_PT_COORD(Crv -> PType);
  165.     CagdRType
  166.     t1 = 1.0 - t;
  167.     CagdCrvStruct *RCrv;
  168.  
  169.     /* Can not compute for linear curves. */
  170.     if (k <= 2) return NULL;
  171.  
  172.     if (APX_EQ(t, 0.0)) {
  173.     /* Use Crv starting tangent direction. */
  174.     CagdCoerceToE3(P1.Vec, Crv -> Points, 0, Crv -> PType);
  175.     CagdCoerceToE3(P2.Vec, Crv -> Points, 1, Crv -> PType);
  176.     CagdCoerceToE3(P3.Vec, Crv -> Points, 2, Crv -> PType);
  177.     }
  178.     else if (APX_EQ(t, 1.0)) {
  179.     /* Use Crv ending tangent direction. */
  180.     CagdCoerceToE3(P1.Vec, Crv -> Points, k - 3, Crv -> PType);
  181.     CagdCoerceToE3(P2.Vec, Crv -> Points, k - 2, Crv -> PType);
  182.     CagdCoerceToE3(P3.Vec, Crv -> Points, k - 1, Crv -> PType);
  183.     }
  184.     else {
  185.     RCrv = BzrCrvNew(k, Crv -> PType);
  186.  
  187.     /* Copy Crv into RCrv, so we can apply the recursive algo. to it. */
  188.     for (i = 0; i < k; i++)
  189.         for (j = IsNotRational; j <= MaxCoord; j++)
  190.         RCrv -> Points[j][i] = Crv -> Points[j][i];
  191.  
  192.     /* Apply the recursive algorithm to RCrv. */
  193.     for (i = 1; i < k; i++)
  194.         for (l = 0; l < k - i; l++)
  195.         for (j = IsNotRational; j <= MaxCoord; j++)
  196.             RCrv -> Points[j][l] = RCrv -> Points[j][l] * t1 +
  197.                        RCrv -> Points[j][l+1] * t;
  198.  
  199.     CagdCoerceToE3(P1.Vec, RCrv -> Points, 0, RCrv -> PType);
  200.     CagdCoerceToE3(P2.Vec, RCrv -> Points, 1, RCrv -> PType);
  201.     CagdCoerceToE3(P3.Vec, RCrv -> Points, 2, RCrv -> PType);
  202.  
  203.     CagdCrvFree(RCrv);
  204.     }
  205.  
  206.     CAGD_SUB_VECTOR(P1, P2);
  207.     CAGD_SUB_VECTOR(P2, P3);
  208.  
  209.     CROSS_PROD(P3.Vec, P1.Vec, P2.Vec);
  210.  
  211.     if ((t = CAGD_LEN_VECTOR(P3)) < EPSILON)
  212.     return NULL;
  213.     else
  214.     CAGD_DIV_VECTOR(P3, t);                /* Normalize the vector. */
  215.  
  216.     return &P3;
  217. }
  218.  
  219. /******************************************************************************
  220. * Return a normalized vector, equal to the normal to Crv at parameter t.      *
  221. * Algorithm: returns the cross product of the curve tangent and binormal.     *
  222. ******************************************************************************/
  223. CagdVecStruct *BzrCrvNormal(CagdCrvStruct *Crv, CagdRType t)
  224. {
  225.     static CagdVecStruct N, *T, *B;
  226.  
  227.     T = BzrCrvTangent(Crv, t);
  228.     B = BzrCrvBiNormal(Crv, t);
  229.  
  230.     if (T == NULL || B == NULL) return NULL;
  231.  
  232.     CROSS_PROD(N.Vec, T -> Vec, B -> Vec);
  233.  
  234.     CAGD_NORMALIZE_VECTOR(N);                /* Normalize the vector. */
  235.  
  236.     return &N;
  237. }
  238.  
  239. /******************************************************************************
  240. * Return a new curve, equal to the derived curve. If the original curve is    *
  241. * rational, NULL is return (curve must be non-rational for this one).          *
  242. * Let old control polygon be P(i), i = 0 to k-1, and Q(i) be new one then:    *
  243. * Q(i) = (k - 1) * P(i+1) - P(i), i = 0 to k-2.                      *
  244. ******************************************************************************/
  245. CagdCrvStruct *BzrCrvDerive(CagdCrvStruct *Crv)
  246. {
  247.     CagdBType
  248.     IsNotRational = !CAGD_IS_RATIONAL_CRV(Crv);
  249.     int i, j,
  250.     k = Crv -> Length,
  251.     MaxCoord = CAGD_NUM_OF_PT_COORD(Crv -> PType);
  252.     CagdCrvStruct *DerivedCrv;
  253.  
  254.     if (!IsNotRational || k < 3)
  255.     FATAL_ERROR(CAGD_ERR_RAT_LIN_NO_SUPPORT);
  256.  
  257.     DerivedCrv = BzrCrvNew(k - 1, Crv -> PType);
  258.  
  259.     for (i = 0; i < k - 1; i++)
  260.     for (j = IsNotRational; j <= MaxCoord; j++)
  261.         DerivedCrv -> Points[j][i] =
  262.         (k - 1) * (Crv -> Points[j][i+1] - Crv -> Points[j][i]);
  263.  
  264.     return DerivedCrv;
  265. }
  266.  
  267. /******************************************************************************
  268. * Convert a Bezier curve into Bspline curve by adding an open knot vector.    *
  269. ******************************************************************************/
  270. CagdCrvStruct *CnvrtBezier2BsplineCrv(CagdCrvStruct *Crv)
  271. {
  272.     CagdCrvStruct *BspCrv = CagdCrvCopy(Crv);
  273.  
  274.     BspCrv -> Order = BspCrv -> Length;
  275.     BspCrv -> KnotVector = BspKnotUniformOpen(BspCrv -> Length,
  276.                                BspCrv -> Order, NULL);
  277.     BspCrv -> GType = CAGD_CBSPLINE_TYPE;
  278.     return BspCrv;
  279. }
  280.